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ABSTRACT 

The  architecture  for  a  new  large-scale  (5  to  100  km,  1  hour  to  1  day)  sea  ice  dynamics 
model  based  on  an  anisotropic  constitutive  law  is  presented  here.  This  architecture  accounts 
directly  for  reffozen  lead  systems  in  the  pack  ice  strength  (with  an  anisotropic  failure  surface)  and 
in  the  ice  thickness  distribution  (with  an  oriented  thickness  distribution).  The  lower  limit  (5  km) 
of  the  model  resolution  is  controlled  by  the  fracture  spacing  of  old,  thicker  ice  and  the  maximum 
lead  width.  The  upper  limit  of  the  model  resolution  (100  km)  is  controlled  by  curvature  in  the 
lead  directions  and  variations  in  the  lead  width.  These  in  turn  are  controlled  by  the  variations  in 
internal  ice  stress  due  to  driving  forces  (winds  and  currents),  which  set  the  time  resolution.  This 
architecture  features  abrupt  changes  in  the  failure  surface  and  the  associated  flow  rule  that  cannot 
be  averaged  over  a  time  step.  In  addition,  the  principal  stress  normal  to  a  new  lead  must  be  zero 
as  it  opens. 

This  model  has  sub-scale  simulations  that  allow  for  the  inclusion  of  phenomena  such  as 
ridging,  rafting,  buckling,  and  fracture  on  the  behavior  of  the  ice.  With  this  new  ice  constitutive 
law,  it  is  possible  to  directly  test  the  ice  failure  strength,  plastic  flow  rule,  and  ice  thickness 
distribution.  The  data  most  useful  for  this  testing  come  from  ice  stress  and  position  buoys 
together  with  SAR  deformation  data.  Some  data  comparisons  have  already  been  made. 

INTRODUCTION 

The  ice  model  constitutive  relations  developed  in  this  paper  describe  the  anisotropic 
behavior  of  pack  ice  that  is  exhibited  by  the  formation  and  evolution  of  leads  in  rafts  and  ridges. 
This  model  is  intended  for  use  when  ice  stress  is  important  and  has  a  marked  influence  upon  the 
ice  dynamics,  which  are  the  nine  or  ten  months  of  the  Arctic  winter.  During  the  summer,  there  is 
very  little  pack  ice  stress,  the  dynamics  are  little  affected  by  the  stress,  and  therefore  there  is  no 
need  for  a  constitutive  law.  The  stress  and  deformation  state  of  a  two-dimensional  element  of 
pack  ice  is  described  by  three  stress  components  in  the  plane  ( ,  cjyy ,  )  (two  normal  and 

one  shear)  and  three  deformation  rates  (two  normal  and  one  shear). 

In  early  ice  dynamics  models  (Coon  et  aL,  1974;  Pritchard,  1975;  Thorndike  et  aL,  1975; 
Coon,  1980;  Pritchard,  1981;  Hibler,  1986),  the  behavior  of  the  ice  was  taken  to  be  isotropic. 
Therefore,  the  behavior  was  represented  by  two  stress  invariants  and  two  deformation  invariants. 
In  these  models,  the  ice  is  considered  as  a  plastic  material  that  has  an  elastic  or  viscous  behavior 
for  stresses  below  failure.  The  isotropic  behavior  is  modeled  by  a  yield  function  with  fixed. 
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selected  shape  in  the  two  stress-invariant  space.  The  failure  strength,  Le.,  size  of  the  yield 
function,  is  controlled  by  one  parameter,  its  compressive  strength.  Such  a  model  can  best 
describe  pack  ice  when  there  are  no  active  leads  (no  plastic  deformation)  or  many  (three  or  more) 
active  lead  directions  in  an  element.  Except  near  fixed  boundaries,  element  sizes  of  hundreds  of 
kilometers  and  time  steps  of  weeks  are  required  to  have  three  or  more  active  lead  directions. 
Although  many  leads  may  show  up  in  an  aerial  photo  or  satellite  image,  only  a  few  are  active. 

Figure  1  in  Stem  et  aL  (1995),  a  typical  ice  motion  product  derived  from  SAR  images, 
shows  the  relative  ice  motion  during  a  three-day  period  in  winter  for  an  80x300  km  area.  The 
relative  ice  motion  is  displayed  as  an  array  of  vectors  on  a  5-km  grid  indicating  the  deviations 
from  the  mean  ice  displacement,  or  is  displayed  as  the  deformed  configuration  of  an  initially 
regular  5-km  grid.  There  are  clearly  several  regions  of  uniform  motion  separated  by  sharp 
discontinuities.  Within  the  regions,  most  of  the  cells  are  rigid,  with  little  or  no  deformation  or 
area  change.  The  discontinuities  between  the  regions  are  essentially  linear  for  the  100  kilometers 
seen  in  the  figure.  Therefore,  isotropic  models  are  not  very  good  approximations  of  ice  behavior 
on  scales  of  5  to  100  kilometers  with  time  resolutions  of  fractions  of  a  day  to  a  few  days.  Much 
of  the  important  physics  is  controlled,  however,  by  leads,  which  need  to  be  resolved  on  time 
scales  that  are  fractions  of  a  day. 

Coon  et  aL  (1992)  suggested  an  anisotropic  model  to  overcome  these  shortcomings  of  the 
isotropic  model  In  the  present  model  two  types  of  ice  behavior  are  considered.  In  the  first 
behavior,  the  deformations  occur  as  a  combination  of  elastic  stretching  and  plastic  stretching  of 
the  one  or  two  existing  lead  directions.  In  the  second  ice  behavior,  one  or  two  leads  are  created 
and  the  deformations  occur  as  a  combination  of  elastic  stretching  and  plastic  stretching  of  the  new 
and  existing  leads.  This  model  determines  the  new  lead  directions.  In  both  behaviors,  the  change 
of  ice  thickness  distribution  in  the  lead  directions  is  calculated. 

A  large-scale  ice-dynamics  model  can  be  described  by  considering  a  momentum  balance, 
driving  forces,  and  a  constitutive  law  for  the  ice  (see  Coon  et  al,  1974).  In  the  application 
considered  for  this  model,  resolutions  of  5  to  100  km  and  1  hour  to  1  day  are  required.  The  time 
and  space  scales  that  can  be  resolved  by  this  model  are  controlled  by  the  driving  forces  and  the 
constitutive  equation.  To  accomplish  this  resolution,  the  pack  ice  is  divided  into  two  types,  as 
illustrated  in  Figure  1: 

•  Isotropic  ice  is  consolidated  multi-year  and  heavy  first-year  ice  with  numerous  flaws 
(brine-fiUed  cracks)  which  limit  its  strength.  Isotropic  ice  does  not  consist  of  single  or 
multiple  floes  but  of  large  clusters  (1  to  20  km)  of  small,  interlocking  rigid  fragments 
which  are  only  100  to  500  m  in  extent. 

•  Oriented  ‘ice’  includes  (following  the  nomenclature  of  WMO,  1970)  open  water,  new  ice, 
nilas,  young  ice,  first- year  ice,  rafted  ice,  and  ridged  ice,  when  they  occur  in  long,  narrow 
features.  First-year  ice,  rafted  ice,  and  ridged  ice  are  transferred  to  isotropic  ice  as  soon 
as  their  strengths  exceed  that  of  isotropic  ice.  When  not  appearing  in  long,  linear  features, 
aU  these  ice  types  are  isotropic  ice.  In  this  paper,  the  term  ‘lead’  is  used  broadly  to 
describe  these  long,  narrow  features  whether  or  not  they  contain  open  water,  new  ice, 
nilas,  young  ice,  first-year  ice,  rafted  ice,  or  ridged  ice.  There  may  or  may  not  be  an  ice 
velocity  discontinuity  across  the  lead.  Furthermore,  parallel  leads  (lead  systems)  are  not 
distinguished  from  a  single  lead  in  this  model 

Table  1  lays  out  the  architecture  of  this  anisotropic  elastic-plastic  sea  ice  mechanics  model. 
The  three  sections  describe  aspects  of  the  model:  ice  dynamics  model,  isotropic  ice  behavior,  and 
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oriented  ice  behavior.  Each  aspect  includes  a  column  of  the  type  of  data  used  to  study  each 
aspect,  a  list  of  the  features  of  each  aspect,  and  the  sub-scale  models  needed  to  describe  each 
aspect.  The  ice  dynamics  aspect  is  essentially  similar  to  previous  large-scale  sea-ice-dynamics 
models.  The  isotropic  ice  behavior  is  substantially  changed;  deformations  involve  leads.  The 
oriented  ice  behavior  is  new  to  pack-ice  constitutive  laws.  Ice  velocity  solutions  for  the 
anisotropic  problem  are  obtained  in  exactly  the  same  way  as  they  are  for  isotropic  problems 
except  that  the  constitutive  law  is  different.  The  interpretation  of  the  ice  velocity  field  at  each 
time  step  is  now  made  in  terms  of  leads  to  generate  the  new  yield  surface  for  the  next  time  step. 
This  procedure  is  very  similar  to  changing  the  size  of  the  yield  surface  in  an  isotropic  model  after 
each  time  step. 

In  the  next  section,  the  behavior  of  isotropic  ice  is  examined,  which  leads  into  a  discussion 
on  the  formation  of  new  leads.  A  section  on  the  behavior  of  oriented  ice  follows,  including  how 
the  anisotropic  yield  surface  is  formed,  the  corresponding  kinematics,  and  how  the  oriented  ice 
thickness  is  redistributed  due  to  deformations  (which  is  the  topic  of  a  companion  paper,  Coon  et 
aL,  1997a). 

ISOTROPIC  ICE 

•  Isotropic  ice  has  an  isotropic  failure  surface  that  remains  relatively  constant  throughout 

the  winter  season.  This  failure  surface  can  be  written  in  terms  of  principal  stresses. 

Isotropic  ice  is  the  ice  surrounding  active  lead  systems;  therefore,  it  may  have  a  wide 
distribution  of  ice  thickness  and  includes  all  multiyear  ice  and  heavy  first-year  ice.  Evans  and 
Untersteiner  (1971)  showed  how  thermally  induced  ice  stresses  will  cause  bending  failure  and 
vertical  cracks  spaced  at  100  to  300  meters,  depending  upon  the  ice  thickness.  As  illustrated  in 
Figure  2,  isotropic  ice  is  made  up  of  large  clusters  of  small,  interlocking,  closely  fitting,  angular, 
rigid  fragments  which  are  only  100  to  500  m  in  size  because  some  of  the  thermal  cracks  flood  and 
freeze.  The  cracks  between  the  fragments  are  not  flooded  with  sea  water  and  offer  a  small 
amount  of  cohesion  due  to  the  material  remaining  below  the  crack.  The  small-scale  mechanisms 
that  form  the  isotropic  ice  and  control  its  large-scale  properties  are  thus  crack  formation,  brine 
drainage  into  the  cracks,  and  partial  freezing  of  the  brine  if  temperatures  drop. 

Since  the  isotropic  ice  consists  of  numerous  fragments  as  shown  in  Figure  2,  the  new  leads 
are  never  as  straight  as  shown  in  Figure  1.  As  illustrated  in  Figure  2  by  the  bold  lines,  they  will 
progress  around  the  fragments  following  a  general  direction  through  the  cells.  Golombek  and 
Banerdt  (1990)  discuss  the  allowable  deviation  of  an  acceptable  crack  direction  from  the  preferred 
lead  direction.  A  deviation  of  15®  from  the  preferred  direction  requires  only  a  modest  increase  in 
stress  at  failure  (less  than  20%  increase);  the  stress  at  failure  is  still  a  factor  of  ten  less  than  that 
for  a  contiguous  material 

Description  of  isotropic  ice  strength 

The  isotropic  ice  is  a  homogeneous  and  plastic  material  with  in-plane  and  out-of-plane 
failure  mechanisms  analogous  to  the  Coulomb  failure  mechanisms  in  three-dimensional  materials. 
The  Coulomb  failure  mechanism  is  well  characterized  in  the  literature  for  soils  (e.g.,  Bowles, 
1968),  for  sea  ice  at  the  laboratory  scale  (e.g.,  Coon  et  al.,  1984),  and  for  pack  ice  (e.g., 
Pritchard,  1988).  The  Coulomb  failure  criterion  states  that  the  shear  stress  on  any  surface  within 
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the  material  must  not  exceed  the  shear  strength  of  that  surface.  The  shear  strength,  in  turn, 
depends  linearly  upon  the  normal  stress  on  that  surface. 

We  can  evaluate  the  strength  parameters  for  isotropic  ice  using  sea-ice  stress  data  such  as 
that  obtained  during  the  U.S.  Navy  Office  of  Naval  Research  (ONR)  Sea  Ice  Mechanics  Initiative 
(SIMI)  field  program  (see,  for  example,  Coon  et  aL,  1997b).  Figure  3  shows  SIMI  stress  data 
plotted  as  the  shear-resultant  invariant  (half  the  difiference  of  the  principal  stress-resultants) 
plotted  against  the  pressure-resultant  invariant  (average  of  the  principal  stress-resultants).  The 
data  points  in  Hgure  3  are  the  hourly  averaged  sea-ice  stress  measurements  made  by  the  authors 
during  the  SIMI  program  in  1993-94  in  the  Beaufort  Sea,  north  of  Alaska.  We  plotted  the  stress 
data  in  search  of  the  envelope  of  the  ice-stress  data  that  corresponds  to  the  failure  surface  of  the 
pack  ice.  That  is,  pack-ice  stress  can  have  a  wide  range  of  values  within  the  envelope,  but,  in  a 
situation  where  it  contacts  the  envelope,  the  pack  ice  fails.  It  forms  a  slip  line  or  opens  a  lead 
system,  thus  limiting  the  stress. 

The  slopes  and  intercepts  of  the  failure  surface  were  chosen  to  enclose  most  of  the  stress 
data  points.  This  failure  surface  is  well  represented  by  a  Coulomb  failure  surface  isotropic  in 
three  dimensions  with  an  unconfined  compressive  strength  of  250  kN/m  and  an  unconfined  tensile 
strength  of  50  kN/m.  The  line  representing  the  in-plane-shear  failure  is  the  middle  line  segment 
shown  in  Figure  3.  The  corresponding  value  for  the  fiiction  angle  is  42°  and  for  the  cohesion  is 
56  kN/m.  These  coefficients  correspond  to  heavy  pack  ice  with  both  first- year  and  multi-year  ice. 
The  unconfined  compressive  and  tensile  strengths,  which  are  material  constants,  define  this  failure 
criterion.  It  is  assumed  that  the  cracks  between  the  floe  fragments  control  the  strength  of  the 
isotropic  ice. 

To  describe  the  stress  state  of  pack  ice,  we  introduce  the  stress-resultant  as  the  product  of 
the  stress  and  the  thickness  of  the  ice.  This  is  a  plane-stress  problem  with  Gj  and  Oj  being  the 

principal  stress  resultants  in  the  plane  of  the  ice,  and  the  stress  resultant  normal  to  the  ice,  a,  , 

being  zero.  Failure  occurs  when  the  stress-resultant  state  is  on  the  failure  surface.  Figure  4 
shows  an  isotropic  failure  surface  in  principal  stress-resultant  coordinate  system.  The  material 
properties  that  define  this  failure  surface  are  the  cohesion  C  and  fiiction  angle  p.  In-plane  failure 
in  the  horizontal  physical  plane  is  represented  by  the  two  lines  in  the  principal-stress  coordinate 
system  marked  in  Figure  4  as  1  and  2.  The  equations  for  these  two  failure  surfaces  are 

Oj  (1  -t-  sinp)  =  a2  (1  -  sinP)  -i-  2Ccos  P  (1) 

^2  (1  +  sinP)  =  Gj  (1  -  sinP)  +  2Ccos  p  (2) 

In  addition  to  these  two  failure  surfaces,  there  are  four  other  failure  surfaces  representing 
failures  out  of  the  plane  of  the  ice  sheet  marked  in  Figure  4  as  3,  4,  5,  and  6.  These  correspond  to 
tensile  failure  (due  to  out-of-plane  shear)  and  compressive  failure  (due  to  out-of-plane  shear);  the 
latter  creates  pressure  ridges  in  thick  ice.  The  equations  for  these  four  failure  surfaces  are 


Gj  =Gt 

(3) 

G2  = 

(4) 

Gi  =  -Gc 

(5) 
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C2  =  -Gc  (6) 

After  comparing  Rgures  3  and  4,  we  can  suggest  other  failure  surfaces.  The  data  in 
Figure  3  show  an  apparent  maximum  value  for  shear  stress-resultant  invariant  of  about  100  kN/m, 
corresponding  to  a  Tresca  or  von  Mises  failure  surface.  This  feature  is  shown  in  Figure  5a  as  an 
additional  failure  surface  that  is  parallel  to  the  pressure  invariant  axis;  it  corresponds  to  an  in¬ 
plane,  pure-shear  failure  that  is  independent  of  the  normal  stress-resultant  and  that  generates  no 
open  water.  The  prior  failure  surfaces  in  Figure  5  are  shown  as  tensile  failure,  in-plane  shear 
failure,  and  compressive  failures.  Hgure  5b  shows  an  alternate  with  zero  tensile  strength;  moving 
the  tensile  failure  surfaces  to  the  principal  stress-resultant  axes  is  a  reasonable  approximation  of 
the  stress  data  and  will  reduce  difficulties  anticipated  in  numerical  solutions. 

FORMATION  OF  NEW  LEADS 

It  is  difficult  to  describe  detailed  physical  mechanisms  for  the  failure  surfaces  shown  in 
Figures  3-5  and  to  determine  corresponding  flow  rules.  It  is  clear  that  new  leads  result  from  the 
failure  of  isotropic  ice  as  well  as  the  condition  that  open  water  can  support  no  normal  stress  or 
shear  stress.  Therefore,  a  lead  cannot  form  when  there  is  compression  in  both  principal 
directions.  When  one  principal  stress  is  zero,  open  water  can  be  created  in  a  lead  parallel  to  the 
non-zero  principal  stress.  The  failure  surface  configuration  shown  in  Figure  5b  allows  tensile 
failure  to  create  open  water  in  a  lead  parallel  to  the  non-zero  principal  stress-resultant  direction 
with  a  zero  principal  stress-resultant  normal  to  the  lead. 

The  failure  surface  in  principal  stress  coordinates  represents  the  failure  of  isotropic  ice,  but 
it  does  not  address  the  yield  surface  or  the  flow  rule.  The  stress  state  on  a  new  lead  must  be 
examined  to  determine  the  yield  surface.  Figure  6  shows  the  lead  coordinates  as  a  set  of  (^,  t|) 
axes  rotated  6  radians  from  the  reference  directions  (x,  y).  There  is  no  traction  across  an  open 
lead,  so  therefore  the  principal  stresses  are  aligned  with  the  lead  and  the  principal  stress  normal  to 
the  lead  is  zero.  For  example,  if  we  have  ai  =  0,  then  the  stresses  are  known  to  be 

=  <^2  =  0  <7^  =  0  (7) 

The  stress  and  deformation  of  a  lead  at  its  creation  are  reported  by  Coon  et  aL  (1995a). 
Simultaneous  measurements  of  sea  ice  stress  on  both  sides  of  a  lead,  ice  motion  from  sequential 
SAR  imagery,  and  ice  motion  from  drifting  GPS/Argos  buoys  are  compared.  The  lead,  which 
appeared  early  in  Day  35  (Feb.  4),  1994,  was  observed,  relative  motions  were  calculated  across 
and  along  the  lead  direction,  and  stresses  were  transformed  to  this  same  coordinate  system.  The 
lead  orientation  and  approximate  location  were  determined  from  SAR  images  of  the  pack  ice 
surrounding  the  SIMI  camp,  as  well  as  net  ice  motion  at  3-day  intervals.  The  ice  motion  around 
the  lead  was  provided  in  greater  temporal  detail  by  an  array  of  drifting  GPS/Argos  buoys 
producing  hourly  positions.  Sea  ice  stress  was  measured  on  each  side  of  the  lead  using  flat-jack 
sensors;  Figure  7  shows  the  stress  components  in  lead  coordinates. 

Shortly  before  the  lead-opening  motion  is  apparent  in  the  ice  motion  data,  the  ice-stress 
state  abruptly  changes  from  bi-axial  compression  to  uni-axial  compression  with  traction  on  the 
lead  near  zero.  The  stress  path  in  stress-resultant-invariant  space  is  shown  in  Figure  8.  After  the 
lead  is  opened,  the  principal  directions  of  strain  rate  and  stress  no  longer  coincide  (Figure  9). 
Therefore,  pack  ice  with  an  active  lead  is  a  highly  anisotropic  material;  isotropic  models  cannot 
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describe  this  behavior  since  they  assume  that  the  principal  directions  of  stress  and  stretching 
always  coincide. 

In  stress  space,  the  yield  is  occurring  on  the  Ca  axis.  As  long  as  there  is  open  water  in  the 
new  lead,  the  stress  state  must  be  represented  by  Equation  (7),  and,  therefore,  the  lead  opening 
and  shearing  motions  represent  deformations  on  a  new  yield  surface  defined  by  Equation  (7). 
This  yield  surface  is  anisotropic  and,  in  Crx-<5yy-<5xf  space,  the  equations  for  the  stress  are  a 
parametric  representation  of  a  line  in  stress  space  in  terms  of  cr^; 

0 

<^xx  =  -|^(i  -  cos2e) 

Gyy  =  -|^(1  +  COS20)  (8) 

rixy  =  — ^(sin26) 

For  any  lead  orientation.  Equation  (8)  is  a  family  of  lines  that  pass  through  the  origin  of 
stress  space,  forming  the  tensile  cutoff  cone  shown  in  Figure  5b. 

ORIENTED  ICE 

•  Oriented  ‘ice’  provides  orientation  to  the  strength  of  the  pack  ice.  Each  lead  has  its  own 

orientation. 

The  oriented  ice  in  leads  and  refrozen  leads  causes  anisotropy  in  the  strength  of  the  pack. 
The  constituents  of  oriented  ice  include  (following  the  nomenclature  of  WMO,  [1970])  open 
water,  new  ice,  nilas,  young  ice,  first-year  ice,  rafted  ice,  and  ridged  ice  when  they  occur  in  long, 
narrow  features. 

Yield  surface  for  oriented  ice  and  flow  rule 

The  generic  shape  of  ±e  lead  yield  surface  is  shown  in  Figure  10a  as  a  rectangular  prism 
with  rounded  ends.  The  prism  represents  the  failure  and  yield  surface  for  the  oriented  ice  in  one 
lead.  The  stress  parallel  to  the  lead  (a^j)  is  aligned  with  the  long-axis  of  the  prism,  or^  is  the  stress 
normal  to  the  lead,  and  is  the  shear  stress  on  the  lead.  When  one  lead  direction  is  active,  the 
anisotropic  yield  surface  for  the  pack  ice  consists  of  the  prism  and  the  rounded  ends  that  are 
contributed  by  the  failure  surface  for  the  isotropic  ice.  When  the  stress  state  lies  on  one  of  these 
ends,  the  isotropic  ice  will  fail  When  the  stress  state  lies  on  the  sides  of  the  prism,  the  lead  fails 
by  a  combination  of  shearing  and  opening  or  closing.  The  top  and  bottom  surfaces  in  Figure  10a 
correspond  to  shear  yield  of  the  oriented  ice  in  the  lead;  the  left  surface  corresponds  to 
compressive  failure  (rafting  and  ridging);  and  the  right  surface  corresponds  to  tensile  failure 
(opening). 

In  this  stress  space,  the  oriented  ice  yield  surface  is  a  prism  of  constant  cross-section  that 
lies  at  45®  to  the  pressure  invariant  axis.  The  pressure  invariant  axis  is  shown  bisecting  the  angle 
between  the  parallel  stress  (a^)  and  the  normal  stress  (a^).  As  an  example,  the  yield  surface  for  a 
lead  with  open  water  is  a  line  lying  along  the  a^-axis,  as  shown  in  Figure  10b.  As  ice  grows  in 
the  lead,  the  cross-section  of  the  yield  surface  grows  as  well  Practically  speaking,  the  precise 
shape  of  that  cross-section  is  less  important  than  its  overall  dimensions  because  the  dimensions 
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change  rapidly  due  to  ice  growth  and  deformation.  Suitable  cross-sections  include  a  rectangle,  a 
diamond,  and  an  ellipse,  as  illustrated  in  Figure  10c.  Future  work  may  determine  whether  it  is  a 
rectangle,  diamond,  oval,  ellipse,  or  other  shape,  but  each  yield  surface  will  have  a  maximum 
stress  perpendicular  to  the  lead  (tensile  strength),  a  minimum  stress  perpendicular  to  the  lead 
(compressive  strength),  and  a  maximum  shear  stress  magnitude  (shear  strength),  associated  with  a 
specific  lead  direction.  A  discussion  of  the  yield  strength  of  oriented  ice  is  now  a  matter  of 
defining  these  strength  parameters. 

It  is  difBcult  to  visualize  the  intersection  of  the  lead  yield  surface  and  the  isotropic  failure 
surface.  Figure  lOd  shows  the  =  0  slice  of  the  isotropic  failure  surface  from  Figure  4  in  the 
same  lead- stress  coordinates  used  in  Figure  10a.  The  slice  of  the  lead’s  rectangular  prism  is 
shown  as  two  lines  parallel  to  the  c^-axis.  When  =  0,  the  other  two  lead-stress  coordinates 
(a^  and  a^)  are  the  principal  stresses  and  the  isotropic  failure  surface  lies  in  the  (a^^-a^)  plane. 

When  there  are  two  oriented  yield  surfaces,  then  the  isotropic  failure  surface  may  not 
come  into  play  at  aU.  For  example,  if  there  two  leads  with  open  water  that  are  not  parallel,  then 
each  contributes  a  line  yield  surface  and  their  intersection  at  the  origin  is  the  yield  surface, 
implying  a  null  stress  state  and  free  drift  of  the  pack  ice.  A  lead  with  open  water  combined  with  a 
lead  with  thin  ice  will  generate  a  yield  surface  consisting  of  a  line  chopped  off  by  the  rectangular 
prism  of  the  other  lead.  As  ice  grows  in  both  leads,  both  yield  surfaces  grow  in  size.  The 
rectangular  prism  grows  beyond  the  isotropic  failure  surface,  and  the  line  grows  to  a  rectangular 
prism.  Perhaps  much  of  the  time,  the  ridging  strengths  of  several  old  leads  limit  the  compressive 
strength  of  the  pack  ice,  not  the  isotropic  compressive  strength  discussed  in  Section  2. 

Failure  mechanisms  of  oriented  ice 

Consider  the  rectangular  lead  yield  surface  in  Figure  10c.  As  shown  in  Figure  11,  each 
side  of  the  rectangle  represents  a  lead-ice  failure  mechanism:  opening,  closing,  the  opposite  side 
of  the  lead  moving  to  the  right  (positive  shear),  or  the  opposite  side  moving  to  the  left  (negative 
shear).  The  dimensions  of  the  rectangle  depend  only  on  the  thickness  of  the  thinnest  ice  within 
the  refrozen  lead  system.  When  the  stress  resultant  touches  this  rectangle,  plastic  flow  occurs 
within  the  lead  system;  Le.,  the  lead  ice  opens,  closes,  or  shears.  A  tensile  strength,  a^,  a 

compressive  strength,  a^,  and  a  shear  strength,  define  the  four  sides  of  the  rectangle.  As  a 

consequence  of  the  normal  flow  rule,  the  side  of  the  rectangle  with  positive  shear  stress  (the  top 
side  as  shown)  represents  stress  states  with  plastic  flow  of  the  lead  ice  in  pure  positive  shear  (no 
opening  or  closing  of  the  lead)  since  the  outward  normal  is  in  the  direction  of  positive  shear 
stress.  The  right  side  is  pure  opening,  the  bottom  is  pure  negative  shear,  and  the  left  is  pure 
compression.  At  the  top  right  comer,  a  range  of  stretching  corresponds  to  the  fan  of  outward 
normals  between  positive  shear  and  opening.  The  other  comers  have  similar  fans  of  outward 
normals.  Thus,  this  rectangular  yield  surface  allows  aU  plastic  stretching  vectors  in  this  plane.  In 
3-D  stress  space,  the  four  sides  of  the  rectangle  represent  planar  faces  that  extend  indefinitely  in 
the  stress  resultant  direction  parallel  to  the  lead.  The  component  of  stretching  parallel  to  the  lead 
is  zero  and  corresponds  to  a  component  out  of  the  plane  shown  in  Figure  11.  If  necessary,  such 
ice  deformations  are  accomplished  by  deforming  another  existing  lead  system  or  by  contacting  the 
failure  surface  for  the  isotropic  ice  and  creating  another  one  or  two  lead  systems. 
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Form  of  yield  function  in  x-y-z  stress  coordinates 

The  four  plastic  yield  functions  that  define  the  four  sides  of  the  rectangular  prism  are  given 
by  the  following  equations: 

(|)jj  =  aj  -axxCOs26j  -CxySinlSj  -Cj^ 

(j)jc  =  “<?!  +<yxx‘'Os29j  +GxySin29j 
(j)jp  =  — <3ijyjSin26j  +  (TjjyCOs29 j  “^js 
(jiju  =  (JxxSin29j  —  dxyCOs29j  ~<7js 

where  9j  is  the  angle  from  x-axis  to  the  direction  parallel  to  the  j-th  lead  and  where  the  subscripts 
t,  c,  p,  and  n  refer  to  the  tension,  compression,  positive  shear,  and  negative  shear  sides  of  the 
rectangle,  respectively.  These  four  sides  form  a  rectangular  prism  that  always  lies  at  a  45°  angle 
to  the  Oraxis  (the  axis  of  the  cone  for  the  isotropic  ice  failure  surface).  The  lead  angle  9j  defines 
the  orientation  of  the  j-th  prism  around  the  (Jpaxis.  In  stress-resultant  space,  29j  is  the  angle 
between  the  CrtT'w  plane  and  the  plane  containing  the  Graxis  and  the  j-th  yield  surface  prism. 


Kinematic  equations  for  pack  ice  with  deformations  along  leads 

Coon  et  aL  (1997a)  described  the  kinematics  and  ice-thickness  redistribution  of  pack  ice 
while  accounting  for  lead  and  ridge  orientation.  The  pack  ice  deformations  are  described  in  terms 
of  velocity  discontinuities  rather  than  invariants,  allowing  direct  interpretation  of  open-water 
production  and  ridging.  The  behavior  requires  an  oriented  ice-thickness  distribution  to  account 
for  the  growth  and  ridging  of  ice  in  refrozen  leads. 

Each  active  lead  direction  contributes  to  the  plastic  deformation  of  the  pack  ice.  The  pack 
ice  motion  is  therefore  the  sum  of  the  contributions  of  lead  directions.  It  is  assumed  that  there  are 
negligible  deformations  of  the  isotropic  ice  between  the  leads  (e.g.,  elastic  deformations).  In  this 
model,  the  velocity  derivative  normal  to  the  lead  directions  is  the  only  part  of  the  kinematics  that 
affects  ice-thickness  redistribution  in  the  leads. 

Ice-thickness  redistribution 

The  distribution  of  ice  thickness  at  each  lead  orientation  can  be  described  much  as  others 
have  described  the  distribution  of  ice  thickness  (e.g.,  Thorndike  et  aL,  1975).  Coon  et  aL  (1997a) 
characterized  these  ice-thickness  distributions  assuming  isotropic  ice  has  one  ice-thickness 
category  and  oriented  ice  has  four  ice-thickness  categories  in  N  lead  systems,  and  developed  the 
mathematical  equations  and  calculation  procedures  for  the  ice-thickness  redistribution.  These 
results  could  be  adapted  to  any  number  of  ice  thickness  categories. 

Once  the  ice  thickness  distribution  in  a  lead  system  is  updated,  the  strength  parameters  for 
that  lead  can  calculated  based  on  the  thickness  of  the  thinnest  ice  thickness  category  available. 
The  tensile  strength  is  near  zero.  The  compressive  strength  can  be  estimated  from  theories, 
numerical  models  (Hopkins,  1994),  laboratory  data,  and  field  data  for  buckling  (Coon  et  al., 
1989a),  rafting  (Pritchard  et  aL,  1995),  and  ridging  (Parmerter  and  Coon,  1972)  of  lead  ice.  The 
shear  strength  can  be  estimated  from  buckling  and  the  shear  strength  of  rubble  ice  (Weiss  et  al., 
1981;  Coon  et  aL,  1995b).  All  vary  with  the  ice  thickness.  Given  these  strengths,  the  yield 
surface  associated  with  each  lead  is  combined  with  the  isotropic  ice  failure  surface  into  a  single 
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failure  surface.  It  is  of  interest  to  note  that  the  oriented  ice  yield  surfaces  will  change  quickly  over 
hours  or  days  while  the  isotropic  failure  surface  may  t5e  essentially  constant  for  a  whole  season. 

SUMMARY  AND  CONCLUSIONS 

The  anisotropic  sea-ice  model  describes  the  behavior  of  individual  lead  systems  within  a 
general  framework  of  plasticity.  The  yield  surface  and  normal  flow  rub  of  each  lead  system  is 
viewed  in  coordinates  aligned  with  it.  This  view  is  two  dimensional  because  only  the  traction 
across  the  lead  affects  its  behavior  until  surrounding  isotropic  ice  fails.  When  a  lead  first  opens, 
its  yield  surface  is  a  line  in  stress  space;  that  is,  it  cannot  sustain  normal  or  shear  stress  across  it. 
Its  thickness  distribution  evolves  thermally.  When  thicker  ice  is  located  in  the  lead,  its 
compressive  and  shear  strengths  are  larger.  These  strengths  can  be  estimated  from  selected 
processes  such  as  buckling,  shearing,  rafting,  or  ridging. 

The  failure  surface  of  the  surrounding  isotropic  ice  is  a  body  of  revolution  in  stress  space, 
or  a  polygon  in  stress  invariant  space,  as  in  isotropic  models.  It  is  bounded  by  tensile  cutoff  (a 
principal  stress  is  zero),  a  maximum  shear  strength,  and  a  cap  on  the  compressive  end. 

A  new  lead  may  form  when  at  least  one  principal  stress  is  zero.  Orientation  of  the  new 
bad  is  determined  by  the  orientation  of  principal  direction  of  the  stress  tensor  at  this  time.  There 
is  no  requirement  that  principal  directions  of  stress  and  stretching  be  aligned.  As  the  lead 
subsequently  deforms,  the  direction  of  principal  stretching  is  limited  by  the  restriction  that 
stretching  along  the  bad  direction  is  zero.  This  condition,  however,  does  not  cause  the  direction 
of  principal  stretching  to  be  aligned  with  the  bad.  Once  the  refrozen  bad  ice  has  significant 
compressive  and  shear  strength,  the  principal  direction  of  stress  can  be  arbitrary.  Leads  therefore 
allow  and  promote  anisotropic  behavior  in  pack  ice. 

We  have  assumed  that  the  isotropic  ice  is  fractured  in  all  directions  and  has  zero  tensile 
strength  in  all  directions.  If  the  sea  ice  is  sheared  under  compressive  stresses,  it  deforms  along 
velocity  discontinuity  lines  (closed  and  shearing  bads).  Under  these  conditions,  the  shear  stress 
invariant  is  at  its  maximum  value  (analogous  to  Tresca  and  von  Mises  failure  criteria),  two 
orthogonal  (closed  and  shearing)  bad  directions  exist,  and  they  bisect  the  principal  directions. 
Such  deformations  do  not  cause  funher  damage  to  the  isotropic  material,  and  future  bad 
formation  is  not  affected  by  these  deformations. 

This  anisotropic  model  provides  an  opportunity  to  describe  the  formation  of  open  water 
more  accurately,  even  to  the  extent  of  defining  its  orientation,  which  can  be  compared  directly 
with  new  bads  observed  in  SAR  imagery.  Although  the  anisotropic  constitutive  law  was 
developed  to  describe  behavior  on  smaller  scales,  say  5  km  resolution,  by  incorporating  the 
physical  behavior  from  even  smaller  scales,  it  also  improves  on  our  ability  to  describe  the  one 
variable  that  is  most  critical  to  climate  dynamics  modeling,  the  amount  of  open  water. 

NOTATION 

C  cohesion  (a  stress  resultant)  in  pack  ice  (N/ra) 

X,  y  reference  coordinate  system  (m) 

P  friction  angle  in  pack  ice  (radians) 

0,  0j  in  physical  space,  angle  from  x-axis  to  (j-th)  bad  direction  (radians) 

^,Ti  bad  coordinate  system  (m) 
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ai,<J2 


*^xx  ’  ^yy  ’  ^xy 
<^xx 

^C  » 

<|)jt  »<!>jc  »4>jp 


principal  stress  resultants  (N/m) 

first  and  second  stress  invariant  resultants  (pressure  and  shear  invariant 
resultant,  respectively)  (N/m) 

components  of  stress  resultant  tensor  in  x-y  coordinates  (N/m) 

xx-deviatoric  stress  resultant  in  x-y  coordinates  (N/m) 
components  of  stress  resultant  tensor  in  lead  coordinates  (N/m) 

unconfined  compressive  and  tensile  strength  resultants  of  the  isotropic  ice, 
respectively  (N/m) 

compressive,  shear,  and  tensile  strength  resultant  of  j-th  lead  system, 
respectively  (N/m) 

plastic  yield  function  for  j-th  lead  direction  (N/m) 

plastic  yield  functions  for  tension,  compression,  positive  shear,  and 

negative  shear  for  j-th  lead  direction  (N/m) 
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Table  1.  Architecture  of  Anisotropic  Model 


Data 

Large-scale  sea  ice  dynamics  model 

Sub-scale  models 

Air  presstire 

Air  drag  coefScient 
Ocean  currents 

Ocean  drag 
coefficients 

Air  temperature 

Ice  thickness 

distribution 

Momenmm  balance 

Driving  forces:  air  stress,  water  stress 

Resulting  farces:  Coriolis  forces,  ice  stress 

Resulting  kinematics:  Ice  velocity,  ice  deformation,  ice 
thickness  redistribution 

Atmospheric  boundary 
layer 

Ocean  boundary  layer 
Ocean  currents 

Ice  stress  and  deformation  are  related  through  the 
Constimtive  Law 

Pack  Ice  Constimtive  Law 

•  For  space  scale  5  - 100  km  and  time  scale  1  hr  - 1  day 

•  Elasto-plastic  or  visco-plastic  formulation 

•  Combines  ijorropic  ice  (multiyear  and  heavy  first- 
year)  and  oriented  ice  (1  or  2  lead  directions  provide 
orientation). 

Data 

Isotropic  ice  behavior 

Sub-scale  models 

Stress 

Acoustic 

Fracture 

Temperature 

Isotropic  Ice: 

MY  and  heavy  FY 

Thermal  stress 

Bending  fracture 

Crack  production 

Isotropic  ice  bdiavior: 

Active  cracks  from  thermal  causes  -300  meter 
fragments. 

In-plane  and  out-of-plane  shear  failure  mechanisms 

Lead  orientation 

Data 

Oriented  ice  behavior 

i 

Sub-scale  models 

Stress 

Ridge  profiles 

Ice  ftiction 

Lead  deformation 
Brine  rejection 

Oriented  ice: 

Isotropic  ice  with  1  or  2  lead  directions 

Ice  buckling,  rafting, 
ridging 

Noise  generation 

Oriented  ice  behavior: 

Ice  growth  in  leads 

Yield  surface  for  leads  and  ridges 

Ice  redistribution 
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Figure  1.  Region  of  pack  ice  with  a  wide  refrozen  lead  and  a  narrow  new  lead  forming 
in  the  isotropic  ice. 
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Pressure  resultant  invariant  (kN/m) 


Figure  3.  Coulomb  failure  surface  deteraiined  from  upper  bound  of  measured  ice  stress-resultant  (from 

Coon  ei  al.,  1997b) 


Tensile  failure 


(a)  With  isotropic  pure-shear  failure  surface  added. 


(b)  With  tensile  strength  reduced  to  zero  and  isotropic  pure-shear  failure  surface  added  to  Figure  4. 
Figure  5.  Alternate  isotropic  failure  surfaces  for  pack  ice 


x,y 

U,  V 
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AU,  AV 


Reference  coordinate  system  (m) 

Velocity  components  in  reference  coordinate  system  (m/s) 
Cell  dimensions  in  reference  coordinate  system  (m) 

Angle  from  x-axis  to  j-th  lead  direction  (radians) 

Lead  coordinate  system  (m) 

Velocity  components  in  lead  coordinate  system  (m/s) 

Cell  dimensions  in  lead  coordinate  system  (m) 

Velocity  jumps  across  At^  in  lead  coordinate  system  (m/s) 


Figure  6.  Coordinate  system  and  velocity  components  of  a  rotated  cell  aligned  with  the  lead. 


Shear  stress  on  lead 
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Figure  7.  Stress  components  at  Buoy  2  on  SIMI  Floe  4  expressed  in  lead  coordinates. 
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Figure  8.  Stress  invariants  while  lead  was  active  (Days  35-40) 
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Figure  9.  Comparison  of  principal  directions  for  stress  and  stretching 
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(c)  Suitable  cross-sections  for  the  lead  yield  surface  may  be  a  rectangle  (used  here),  a  diamond, 
an  ellipse,  etc. 


Slice  of  lead  yield 


(d)  Slice  of  isotropic  failure  surface  and  lead  yield  surface  at  zero  lead  shear  stress 


Figure  10.  Oriented  ice  failure  surfaces  (continued) 
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Figure  11.  Oriented  ice  failure  surfaces  with  associated  normal  flow  rule 
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ABSTRACT 

This  paper  describes  the  determination  of  geophysical  pack-ice  stress  from  sea-ice  stress 
measurements  made  during  the  U.S.  Navy  Office  of  Naval  Research  Sea  Ice  Mechanics  Initiative 
field  program  and  describes  the  interpretation  of  the  data  in  terms  of  a  pack-ice  failure  surface. 
This  effort  involved  two  steps:  learning  how  to  interpret  sensor  oil  pressure  as  ice  stress,  and 
learning  under  what  conditions  the  measured  ice  stress  represents  the  geophysical  stress.  These 
data  are  shown  consistent  with  the  indentation  pressure-area  curve,  which  represents  the  strength 
of  ice  over  scales  from  laboratory  samples  to  pack  ice  elements  of  hundreds  of  kilometers. 

INTRODUCTION 

This  paper  reports  on  the  measurement  of  sea  ice  stress  for  many  months  during  the  U.S. 
Navy’s  Office  of  Naval  Research  (ONR)  Sea  Ice  Mechanics  Initiative  (SIMI)  field  program  during 
the  winter  of  1993-1994.  The  autumn- 1994  segment  of  the  SIMI  field  program  is  described  in 
Coon  et  aL  (1994).  Coon  et  al  (1989)  describe  the  stress  sensors  and  the  installation  of  the 
sensors  in  rosettes  to  measure  the  stress  state  in  the  ice.  We  attempted  to  not  only  measure  ice 
stress  but  also  geophysical  pack-ice  stress  resultant.  The  distinction  is  that  the  former  is  a  local 
measurement  (which  varies  with  position  within  a  floe)  and  the  latter  is  an  integral  through  the  ice 
thickness  and  an  average  over  a  few  kilometers.  Geophysical  ice  stress  resultants  are  needed  in 
large-scale  ice-dynamics  models. 

JUSTIFICATION  FOR  INTERPRETING  THIS  DATA  AS  LARGE-SCALE  STRESS 

It  has  been  suggested  that  it  is  not  possible  to  measure  geophysical  pack-ice  stress 
resultant  by  making  a  few  stress  measurements  in  one  floe.  We  present  data  showing  that  stress 
measurements  represent  the  sum  of  local  stresses  and  geophysical  stresses  and  that  geophysical 
stress  can  be  determined.  Our  goals  were  that  (1)  the  oil  pressure  in  the  stress  sensor  could  be 
interpreted  in  terms  of  ice  stress  and  (2)  the  local  stresses  could  be  calculated  and  subtracted  from 
the  measured  stress  to  acquire  the  geophysical  stress.  Three  examples  are  shown  to  support  the 
hypothesis  that  the  sensor  oil  pressures  can  be  related  to  geophysical  ice  stress. 

The  first  example  was  a  ridge-building  event  observed  by  the  scientists  and  crew  aboard 
the  Norwegian  research  ship  Polarbjdm  on  October  7,  1988  (Coon  et  aL,  1989).  At  the  time,  a 
sea-ice  stress  rosette  was  installed  about  300  meters  away.  The  stress  data  show  that  the 
direction  of  the  largest  compressive  principal  stress  points  in  the  direction  perpendicular  to  the 
ridge  when  it  is  building.  Prior  to  the  riding  event,  the  principal  directions  had  large  variations  as 
shown  in  the  data  from  Coon  et  aL  (1989). 

The  second  example  is  the  stress  and  deformation  of  a  lead  at  its  creation  reported  by 
Coon  et  al.  (1995b).  Simultaneous  measurements  of  sea  ice  stress  on  both  sides  of  a  lead,  ice 
motion  from  sequential  SAR  imagery,  and  ice  motion  from  drifting  GPS/ Argos  buoys  are 
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compared.  The  lead,  which  appeared  early  in  Day  35  (Feb.  4),  1994,  was  observed,  relative 
motions  were  calculated  across  and  along  the  lead  direction,  and  stresses  were  transformed  to  this 
same  coordinate  system.  The  lead  orientation  and  approximate  location  were  determined  from 
SAR  images  of  the  pack  ice  surrounding  the  SIMI  camp,  as  well  as  net  ice  motion  at  3-day 
intervals.  The  ice  motion  around  the  lead  was  provided  in  greater  temporal  detail  by  an  array  of 
drifting  GPS/ Argos  buoys  producing  hourly  positions.  Sea-ice  stress  was  measured  on  each  side 
of  the  lead  using  flat-jack  sensors. 

Since  the  calculation  of  sea-ice  stress  from  sensor  oil  pressure  will  affect  all  sensors  in  a 
rosette  equally,  we  calculated  “stress”  components  directly  from  the  sensor  oil  pressures  and  the 
orientation  of  the  sensors.  Figure  1  shows  these  “stress”  components  in  lead  coordinates.  Shortly 
before  a  lead-opening  motion  is  apparent  m  the  ice  motion  data,  the  ice-stress  state  abruptly 
changed  from  bi-axial  compression  to  uni-axial  compression,  with  traction  on  the  lead  near  zero. 
The  stress  path  in  stress  invariant  space  is  shown  in  Figure  2.  The  principal  stress  direction  turned 
until  it  was  nearly  parallel  to  the  new  lead  and  remained  that  way  for  over  a  week.  The  principal 
stress  normal  to  the  lead  dropped  to  zero  just  before  the  lead  opened  and  the  stress  remained  that 
way  for  over  a  week.  The  principal  stress  parallel  to  the  lead  was  variable  but  remained  large. 

The  third  example  is  related  to  determining  the  failure  surfaces  for  pack  ice.  Stress 
analyses  in  Croasdale  et  aL  (1986)  indicated  that  the  stresses  in  the  central  portion  of  the  floe 
reflect  the  average  stress  on  the  floe  perimeter  with  a  reasonable  degree  of  accuracy  for  a  wide 
range  of  loading  conditions.  Hence,  Croasdale  et  aL  suggested  deploying  stress  rosettes  near  the 
center  of  a  large  ice  floe  in  the  pack  ice.  In  essence,  a  large  floe  may  be  used  as  a  transducer  for 
measuring  large-scale  pack  ice  driving  forces.  Coon  et  aL  (1995a)  report  sea-ice  stresses 
measured  in  numerous  locations  m  “Lake  Andy,”  a  large  area  of  flat,  thin  ice  contained  within  the 
very  uneven,  3-km  floe  with  the  SIMI  main  camp  in  the  fall  of  1993.  Even  as  late  as  March  1994, 
Lake  Andy  was  stiU  over  a  kilometer  from  the  edge  of  the  floe.  NWRA  installed  eight  stress 
sensors  in  Lake  Andy  for  a  winter-over  experiment.  AU  of  the  stress  sensors  were  far  from  the 
“shore,”  being  many  characteristics  lengths  for  edge  bending  effects  from  shore  (Hetenyi,  1946). 
These  sensors  should  have  been  exposed  to  the  same  ice  stress  conditions,  whether  they  were 
geophysical  or  thermal  stresses.  Figure  3  (not  to  scale)  shows  the  relative  orientation  of  these 
stress  sensors.  Various  combinations  were  used  to  form  stress  rosettes  of  various  sizes  ranging 
from  3  to  100  meters.  All  of  these  sensors  were  installed  at  35  cm  except  stress  sensor  e,  which 
was  taken  to  be  the  average  of  sensors  at  20  cm  and  50  cm.  Stress  sensors  c,  e,  f,  and  g  are  all 
the  same  orientation.  Along  with  stress  sensors  a  and  b,  these  four  sensors  formed  four  rosette 
combinations  with  diameters  of  3,  10,  50,  and  100  meters.  Again,  we  calculated  stress 
components  directly  from  the  sensor  oil  pressures  and  the  orientation  of  the  sensors;  these 
“stress”  components  are  combined  to  form  the  first  and  second  stress  invariants,  the  pressure  and 
shear  invariants.  A  calibration  factor  affects  all  sensors  equally  and,  therefore,  appears  as  a 
scaling  factor  on  both  invariants.  Thermal  stress  affects  only  the  pressure  invariant  but  should 
affect  aU  sensors  equally  as  weU.  Figures  4-7  show  the  resulting  four  stress-invariant  charts,  using 
stress  invariants  calculated  directly  from  three-hour  averages  of  the  sensor  oil  pressures. 
Although  there  is  some  variation  among  the  four  combinations,  these  four  rosettes  have 
essentially  the  same  stress  envelopes,  making  the  data  useful  for  pack-ice  modeling. 
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SEA-ICE  STRESS  DATA  PROCESSING  PROCEDURE 

NWRA  has  developed  a  sea-ice-stress  data-processing  procedure  that  combines  numerous 
features.  For  the  SIMI  data,  we  have  stress  data  at  one  depth,  temperature  data  at  multiple 
depths,  ice  thickness  sampled  at  a  few  times,  salinity  at  one  time,  and  in-situ  calibration  test  data 
on  many  sensors  but  not  alL  The  sea  ice  data  processing  procedure  developed  by  NWRA  is 
summarized  below. 

1.  Process  in-situ  calibration  test  data  on  the  stress  sensors  taken  at  several  load  rates. 

2.  Interpolate  the  ice  thickness  data. 

3.  Calculate  the  distribution  of  the  elastic  modulus,  E,  through  the  thickness  using  the  ice 
temperature,  salinity,  and  interpolated  ice  thickness  data. 

4.  Remove  the  slow  changes  in  stress  sensor  oil  pressure  resulting  from  the  freeze-in 
stresses,  and  the  sensor  pressurization  cycles. 

5.  Remove  the  long-term  trends  from  ice  temperature  data  to  account  for  creep  in 
thermal  stress. 

6.  Compare  the  de-trended  stress  sensor  oil  pressure  and  the  de-trended  temperature  data 
at  the  sensor  depth  and  remove  sensor  temperature  effects,  if  any. 

7.  Process  the  resulting  stress  sensor  oil  pressure  to  produce  ice  stress  based  on  the  time- 
dependent  model  developed  from  the  in-situ  calibration  tests. 

8.  Calculate  the  thermal  stress  at  stress  sensor  depth  using  the  de-trended  temperature 
data  at  aU  depths  and  adjust  the  magnitude  of  these  thermal  stress  histories  to  fit  the 
measured  stress  during  pure  thermal  events. 

9.  Subtract  the  calculated  thermal  stress  from  the  ice  stress. 

10.  Subtract  the  stress  at  the  depth  of  the  stress  sensor  due  to  the  bending  moment  caused 
by  the  large-scale  stress  and  the  variation  in  the  elastic  modulus  with  depth,  if  any. 

1 1.  Keep  only  the  data  when  at  least  three  sensors  in  a  rosette  have  oil  pressure  over  20 
kPa  for  the  entire  creep-integration  period. 

The  appendix  provides  details  on  the  application  of  this  procedure  to  the  SIMI  stress  data. 

FAILURE  SURFACE  FOR  LARGE-SCALE  STRESS 

Based  on  the  events  and  data  discussed  above,  we  conclude  that  sea-ice  stress 
measurements  represent  large-scale  pack-ice  stress  averaged  over  the  dimension  of  the  typical  floe 
used  to  measure  the  pack-ice  stress,  1  to  3  km.  The  strength  parameters  defining  a  failure  surface 
using  this  stress  data  also  represent  averages  over  1  to  3  km. 

The  outer  envelope  of  the  ice-stress  data  actually  provides  an  upper  bound  for  the  pack- 
ice  failure  surface  because  the  pack  ice  strength  is  often  smaller  than  the  outer  envelope  due  to  the 
presence  of  thin  ice.  This  envelope  represents  the  failure  surface  for  isotropic  pack  ice  as  defined 
by  Coon  et  aL  (1992  and  1997). 

Figure  8  shows  the  SIMI  ice  stress  data  (0.7-m  to  1.6-m  ice)  with  no  correction  for 
thermal  stress.  The  data  shown  on  the  figure  are  hourly  averages  for  only  those  hours  that  have 
three  or  four  sensors  reading  over  20-kPa  oil  pressure  for  the  entire  hour.  This  was  done  because 
our  calibration  tests  and  other  analyses  showed  these  sensors  to  be  an  unreliable  method  to 
measure  tension  and  low  compressive  stresses. 

Figure  9  shows  the  same  SIMI  data  but  with  the  correction  for  thermal  stress.  Since  the 
shear  invariant  is  unaffected  by  the  thermal  stresses,  the  individual  data  points  have  moved 
04/30/97 


4 


horizontally  as  compared  to  Figure  8.  It  is  worthwhile  to  note  that  the  measurements  were  made 
in  a  wide  range  of  ice  thickness:  0.7  to  1.6  m.  These  results  appear  to  be  independent  of  the  ice 
sheet  thickness  in  which  the  sensors  were  installed. 

After  inspection  of  this  sea-ice-stress  data,  we  concluded  that  the  isotropic  failure  surface 
of  the  pack  ice  is  best  characterized  by  a  Coulomb  failure  surface.  It  is  represented  by  three  lines 
on  the  pressure-invariant  plot:  in-plane-shear  failure,  out-of-plane  compressive  failure,  and  out-of¬ 
plane  tensile  failure.  The  Coulomb  failure  mechanism  is  well  characterized  in  the  literature  for 
soils  (e.g.,  Bowles,  1968),  for  sea  ice  at  the  laboratory  scale  (e.g.,  Coon  et  aL,  1984),  and  for 
pack  ice  (e.g.,  Pritchard,  1988). 

Hgure  9  also  shows  lines  indicating  the  recommended  Coulomb  failure  surface.  The 
recommended  failure  surface  is  isotropic  in  three  dimensions  and  has  only  an  unconfined 
compressive  strength  <Jc  and  an  unconfined  tensile  strength  at.  The  corresponding  strength 
parameters  are  a  compressive  strength  of  250  kN/m  and  a  tensile  strength  of  50  kN/m.  These 
strengths  are  consistent  with  those  used  in  large-scale  ice  dynamics  models.  Pritchard  (1993) 
concludes  firom  models  of  observed  events  that  the  unconfined  compressive  strength  must  exceed 
150  kN/m  and  the  shear  strength  must  exceed  27  kN/m. 

Corresponding  values  for  the  Coulomb  cohesion  C  and  the  friction  angle  (])  for  the 
recommended  Coulomb  failure  surface  are  56  kN/m  and  42®,  respectively.  The  steeper  slope  on 
the  right  represents  out-of-plane  tensile  failure;  the  more  gently  sloped  line  on  the  right  represents 
in-plane  shear  failure;  and  the  sloped  line  on  the  left  represents  out-of-plane  compressive  failure. 
The  out-of-plane  compressive  strength  is  the  same  as  the  in-plane  compressive  strength.  This 
failure  surface  represents  the  maximum  pack-ice  strength  at  a  scale  of  1  to  3  km.  It  is  applicable 
to  an  ice  dynamics  model  only  if  the  anisotropic  effects  of  teads  are  incorporated  into  the  model 
This  failure  surface  is  also  applicable  to  determining  ice  loads  on  a  structure  where  it  is  imperative 
to  know  the  maximum  strength  of  ice  on  this  scale. 

COMPARISON  TO  OTHER  DATA 

Figure  10  compares  the  recommended  Coulomb  failure  surface  to  the  largest  yield 
surfaces  used  by  Pritchard  (1980)  and  Hibler  (1980)  in  isotropic  ice-dynamics  models.  The 
isotropic  failure  surfaces  shown  in  the  figure  correspond  to  compact,  heavy  ice  conditions  and 
have  a  maximum  strength  of  100  kN/m.  Under  less  compact  ice  conditions  with  more  thin  ice, 
Pritchard  and  Hibler  use  much  smaller  isotropic  failure  surfaces,  as  discussed  by  Sanderson 
(1988).  Their  isotropic  strengths  were  chosen  such  that  the  velocity  field  in  model  runs  would 
“best”  represent  the  measured  velocity  field  in  some  average  sense.  It  appears  that  the  measured 
stress  data  does  not  fit  inside  the  large-scale  isotropic  yield  surface  very  well. 

The  meaning  of  “large  scale”  is,  however,  very  important  here  and  can  best  be  examined 
by  looking  at  the  strength  of  ice  over  many  scales.  To  this  end,  NWRA  has  developed  a 
procedure  to  use  sea-ice  stress  data  to  add  data  points  to  the  pressure-area  curve  for  ice 
indentation  (Sanderson,  1988).  Figure  9  showed  all  of  our  SIMI  sea-ice  stress  data  plotted  as 
stress-resultant  invariants,  where  the  stress  resultant  is  the  integral  of  the  measured  ice  stress 
though  the  ice  thickness  at  the  stress  measurement  site.  These  ice  stress  data  were  all  measured 
near  the  center  of  ice  floes  that  were  1  to  3  meters  thick  and  1  to  3  km  across.  As  discussed 
above,  these  ice  stress  data  represent  the  average  stress  around  the  perimeter  of  the  instrumented 
floe  and  the  data  shows  a  maximum  compressive  strength  of  250  kN/m. 
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To  compare  this  data  to  the  pressure-area  curve,  we  consider  the  instrumented  floe  to  be  a 
large  indentor.  Coon  et  aL  (1984)  discuss  the  relationship  between  indentation  pressure  and 
material  properties  such  as  compressive  strength  and  friction  angle.  Equation  (17)  in  Coon  et  aL 
(1984)  shows  that  indentation  pressure  equals  the  compressive  strength  for  very  large  values  of 
the  indentor  (floe)  diameter  to  ice  thickness  aspect  ratio.  Since  the  aspect  ratio  in  this  case  is  very 
large,  the  indentation  factor  is  one.  Consequently,  the  maximum  indentation  pressure  is  simply 
the  maximum  compressive  strength. 

These  SIMI  stress  data  can  now  be  plotted  on  Sanderson’s  pressure-area  curve.  For  the 
stated  ranges  in  ice  thickness,  indentor  diameter,  and  maximum  compressive  strength,  the 
indentation  pressures  range  from  0.08  to  0.25  MPa  and  the  contact  areas  range  from  1000  to 
9000  m^.  Figure  1 1  shows  the  excellent  agreement.  Sanderson  has  derived  indentation  pressures 
from  the  isotropic  yield  surfaces  used  by  Pritchard  (1980)  and  Hibler  (1980).  These  are  shown  as 
‘^eso-scale  models”  in  the  lower  right  comer  of  Figure  10  for  maximum  strengths  ranging  from 
5  to  100  kN/m. 

As  a  general  rule,  for  scales  from  0.1  cm  to  500  km,  sea  ice  strength  and  indentation 
pressure  are  controlled  by  three  sea  ice  rheologies  with  three  constitutive  laws: 

1.  For  laboratory  scale  and  intact  sea  ice  (indentor  width  from  0.1  cm  to  10  m),  mechanical 
properties  (e.g.,  elastic  modulus,  compressive  strength)  are  known  to  depend  upon  grain 
size,  brine  volume,  temperature,  and  strain  rates.  Indentation  pressure  depends  upon 
mechanical  properties,  indentor  aspect  ratio,  and  the  specific  failure  mechanism  (e.g.,  in 
plane  or  out-of-plane,  Tresca  or  Coulomb). 

2.  An  ice  floe  is  a  jointed  materiaL  like  igneous  rock.  For  an  indentor  width  from  100  to  500 
m,  the  compressive  strength  is  controlled  by  the  cohesion,  friction,  and  spacing  of  the 
joints  (brine-fiUed  cracks)  occurring  every  100  to  300  m  in  first-year  and  multi-year  sea 
ice.  The  intact  ice  between  joints  does  not  participate  since  the  joints  are  perhaps  100 
times  weaker  than  the  intact  ice. 

3.  Pack  ice  at  the  geophysical  scale  (1  to  500  km),  deformations  occur  along  leads  and  ridges 
spaced  at  10  to  50  km;  the  ice  floes  do  not  participate.  Buckling,  rafting,  and  ridging 
fkilure  mechanisms  (which  depend  upon  the  refrozen  lead  ice  thickness)  control  the 
compressive  strength  and  cause  anisotropy  in  the  pack-ice  strength. 

CONCLUSIONS 

The  SIMI  stress  data  show  that  it  is  possible  to  measure  geophysical  pack-ice  stress.  In 
addition,  we  conclude  that  it  is  inappropriate  to  scale  down  large-scale  isotropic  ice  behavior  to  a 
cell  size  of  a  few  kilometers  without  separately  modeling  the  (anisotropic)  leads  and  ridges. 


APPENDIX  A.  PROCEDURE  FOR  SEA-ICE  STRESS  DATA  PROCESSING 
Calibration  tests 

NWRA  used  thin,  circular,  oil-filled  flat-jack  ice  stress  sensors  to  measure  sea  ice  stress 
during  ONR’s  SIMI  Program  (Coon  et  aL,  1995a;  Coon  et  aL,  1995b).  To  caHbrate  these 
sensors,  we  performed  in-situ  compression  tests  at  two  speeds;  rapid  loading  (around  100 
kPa/rain.),  which  produced  a  square-wave  pressure  history,  and  slow  loading  (around  one 
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kPa/min.),  which  produced  a  ramp  pressure  history.  The  load  rate  of  the  ramp  tests  more  closely 
matched  geophysical  load  rates,  while  the  rapid  load  rate  of  the  square-wave  test  is  fast  enough  to 
measure  the  elastic  response  of  the  ice-sensor  system.  A  large,  flat  air  jack  was  used  to  compress 
the  ice  containing  the  stress  sensor.  The  calibration  tests  were  run  at  separate  sensor  installations 
except  for  one  rosette  that  was  tested  after  a  1-meter  crack  made  the  site  unusable.  Based  on  the 
properties  of  sea  ice  (see  Michel,  1978,  for  example)  and  the  test  results,  we  chose  to  model  the 
response  of  the  sensor  to  ice  stress  with  a  dynamic  inclusion  factor  of  the  form 

<j  =  -AP  -  J  BPdx  (1) 

t-'Co 

where  a  is  the  ice  stress  normal  to  the  sensor  or  the  air  pressure  in  the  air  jack,  P  is  the  oil 

pressure  in  the  sensor,  and  A,  B,  and  Xo  are  constants. 

By  plotting  the  air  pressure  versus  the  oil  pressure  found  in  the  square-wave  tests,  we 
obtained  the  effective  elastic-inclusion  factor.  A,  for  the  sensor  installed  first-year  sea  ice.  In 
many  of  these  plots,  a  change  in  slope  occurs  at  10  to  20  kPa  oil  pressure.  We  judged  that  the 
stress  sensor  readings  are  unreliable  whenever  the  oil  pressure  is  less  than  20  kPa  and,  therefore, 
concluded  that  these  sensors  do  not  measure  tensile  stresses  or  small  compressive  stresses.  The 
creep  terms  in  the  dynamic  inclusion  factor,  B  and  Xo,  were  determined  using  the  linear  up-ramp 
portion  of  each  ramp  test  along  with  the  effective  elastic  inclusion  factor  from  a  prior  square- 
wave  test  at  the  same  installation.  The  measured  air  pressure  for  the  whole  test  was  compared 
with  the  air  pressure  calculated  from  the  oil  pressure  and  the  model 

Although  many  such  square-wave  and  ramp  tests  were  performed,  eight  tests  were 
selected  as  valid  determining  generic  values  for  A,  B,  and  Xo-  Invalid  ramp  tests  were  rejected  for 
any  one  of  the  following  reasons: 

•  The  oil  pressure  never  reached  20  kPa, 

•  The  linear  ramp  time  after  the  oil  reached  20  kPa  was  too  short  (e.g.,  less  than  40 
minutes), 

•  The  oil  pressure  data  was  polluted  by  external  ice  stress  (e.g.  thermal  and  geophysical 
stresses). 

We  found  no  correlation  of  tests  results  to  any  of  the  following  parameters: 

•  The  sensor  diameter  (15  or  20  cm), 

•  Time  of  year  (fall  or  spring  in  first-year  ice), 

•  Air  jack  stand-off  distance  (15  or  30  cm  or  both  on  opposite  sides), 

•  The  measured  elastic  modulus  of  the  stress  sensor. 

The  creep  model  accounts  for  the  effects  of  the  air-jack-loading  rate.  Although  we  could  model 
an  individual  test  weU,  the  standard  deviation  of  A,  B,  and  Xo  for  the  eight  valid  tests  were  37%, 
57%,  and  64%  of  there  values,  respectively.  We  have  chosen  to  use  the  following  values: 

A  =  2.43  B  =  -0.037  min.'*  Xo  =  20  min.  (2) 

Using  these  mean  values,  the  standard  deviation  of  the  calculated  value  of  the  maximum  air 
pressure  during  the  eight  ramp  tests  is  64  kPa. 
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Converting  SIMI  ice  stress  measurements  to  geophysical  stress  resultants 

The  NWRA  SIMI  stress  rosettes  consisted  of  four  sensors  installed  at  45°  compass 
intervals  and  at  roughly  mid-depth  near  the  middle  of  a  large  expanse  of  flat,  first-year  sea  ice. 
Three  rosettes  were  installed  in  the  autumn  in  roughly  70-cm  ice;  one  was  moved  in  the  spring  to 
a  fourth  site  in  1.6-raeter  flat,  first-year  ice.  Time  and  oil-pressure  data  was  stored  on  site  by 
dataloggers  onto  storage  modules  at  5-minute  intervals  and  a  10-minute-interval  data  set  was 
transmitted  to  the  Argos  satellite  system.  Consequently,  the  data  interval  is  five  minutes  for  the 
initial  portion  of  each  data  stream  fi-om  all  four  sites  but  is  ten  minutes  for  the  final  portion  for  the 
three  rosettes  that  were  lost  at  sea. 

To  assure  contact  between  the  sensor  and  the  creeping  sea  ice,  all  sensors  were 
pressurized  at  regular  intervals  by  adding  oil  to  the  sensor  until  the  pressure  increased  by  at  least 
70  kPa.  A  minimum  of  two  hours  of  data  was  lost  for  the  pressurization  cycle.  After 
pressurization,  the  added  pressure  creeps  out,  rapidly  at  first,  then  slower.  We  modeled  this 
process  as  the  sum  of  fast  and  slow  exponential  decays,  fit  each  pressurization  cycle  for  each 
sensor,  and  subtracted  the  double-exponential  curve  firom  the  oil  pressure.  The  two  exponential 
time  constants  were  approximately  three  hours  and  14  days.  The  remaining  oil  pressure  was 
processed  with  the  dynamic  inclusion-factor  model  to  generate  ice  stress.  Whenever  the 
remaining  oil  pressure  was  less  than  20  kPa,  the  data  was  dropped  because  of  our  findings  in  the 
square-wave  calibration  tests  mentioned  above.  Because  Equation  (2)  requires  20  minutes  of 
data,  the  first  20  minutes  after  any  data  dropout  is  lost  as  well 

The  ice  stress  measured  by  the  sensor  is  our  best  estimate  of  the  ice  stress  at  the  depth  of 
the  sensor  in  the  ice  sheet,  but  we  seek  the  geophysical  ice  stress  resultant.  For  that,  we  need  to 
account  for  ice  thickness,  elastic  modulus  variations  through  the  ice  thickness,  thermal  stress,  and 
bending  stress.  The  ice  thickness  determined  firom  drilling  holes  and  fi'om  thermistors  installed 
with  the  stress  rosette  was  consistent  at  the  first  three  sites;  growth  in  ice  thickness  was 
interpolated  in  proportion  to  the  square  root  of  time.  Cox  and  Weeks  (1988)  give  estimates  for 
the  distribution  of  the  elastic  modulus  through  the  ice  thickness  for  various  thicknesses  of  first- 
year  sheet  ice.  We  approximated  their  results  with  a  constant  elastic  modulus  for  the  top  80%  of 
thickness  and  zero  for  the  bottom  20%.  With  this  distribution,  there  is  no  correction  needed  for 
bending  stress  as  long  as  the  sensor  is  in  the  top  80%  of  the  ice  sheet. 

A  theory  for  the  stresses  due  to  temperature  gradients  in  an  unconfined,  thick,  elastic  plate 
of  uniform  thickness  in  which  the  temperatures  and  the  properties  are  a  function  only  of  depth  is 
given  by  Timoshenko  and  Goodier  (1951): 

otEAT  1  12v 

=  <5 w  =  -  -f-  ,  j  aEATdz  — ^  J  oEATzdz  (3) 

“  ^  l-V  h(l-v).|.%  h5(l-v)-i/2 

where  x  and  y  are  the  coordinates  in  the  plane  of  the  plate,  z  is  the  coordinate  down  from  the  mid¬ 
depth  of  the  plate,  a  is  the  stress  at  z,  a  is  the  thermal  expansion  coefficient,  E  is  the  elastic 
modulus,  AT  is  the  temperature  change  from  a  uniform  temperature  as  a  function  of  z,  v  is  the 
Poisson’s  ratio  for  the  material,  and  h  is  the  plate  thickness.  The  depth  used  to  evaluate  the  first 
term  on  the  right  side  of  Equation  (3)  is  the  depth  of  the  resulting  stresses  on  the  left  side.  For 
the  ice  sheet,  we  need  to  evaluate  the  thermal  stress  at  the  depth  of  the  sensor,  a  is  essentially 
constant,  and  E  is  zero  for  the  bottom  20%  of  the  depth.  The  thermal  stress  is  sea  ice  is  thus 
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The  ice  temperatures  were  measured  with  thermistors  at  several  depths.  To  account  for 
the  creep  out  of  the  thermal  stress,  we  calculated  the  temperature  change  not  from  an  initial 
uniform  temperature  but  from  a  reference  temperature  calculated  by  smoothing  the  measured 
temperatures  with  a  low-pass  filter.  The  filter  has  a  14-day  time  constant,  which  is  the  same  time 
constant  as  the  slow  creep  found  in  the  pressurization  decay.  Constant  values  of  1  GPa  for  the 
effective  elastic  modulus,  0.3  for  Poisson’s  ratio,  and  5.1  x  10'^  per  “C  for  the  thermal  expansion 
coefficient  were  found  to  make  the  calculated  thermal  stress,  which  is  the  same  in  aU  (horizontal) 
directions,  consistent  with  the  ice  stress  measured  with  all  four  sensors  in  the  rosettes.  These 
values  are  similar  to  those  used  by  Lewis  (1993)  in  his  calculations  of  thermal  stress  in  sea  ice. 
Lewis,  however,  uses  a  variable  elastic  modulus,  a  different  definition  for  a  reference  temperature, 
and  creep  strain  in  his  model 

The  geophysical  stress  resultant  is  the  integral  of  the  geophysical  stress  through  the 
thickness  of  the  ice  sheet.  Since  the  elastic  modulus  is  constant  for  the  top  80%  of  the  thickness, 
the  difference  of  the  ice  stress  and  the  thermal  stress  is  the  geophysical  stress  for  that  layer.  The 
bottom  20%  must  have  no  stress  since  the  elastic  modulus  is  zero,  so  the  geophysical  stress 
resultant  is  0.8h  times  the  difference  of  ice  and  thermal  stress. 

In  order  to  calculate  the  geophysical  stress-resultant  invariants,  at  least  three  sensors  in  a 
rosette  need  to  have  valid  data  for  20  minutes.  When  that  was  true,  we  calculated  the  invariants 
by  one  of  the  following  two  sets  of  equations.  When  sensors  a,  b,  and  c  were  valid.  Equations  (5) 
and  (6)  provide  the  invariants: 


<^1  =^K+<yc) 


(5) 


(6) 


where  Ca,  <Tb,  and  <3c  represent  the  geophysical  stress  resultant  derived  from  three  sensors  spaced 
at  45°  intervals.  Similar  equations  are  used  for  the  other  sets  of  three,  as  appropriate;  the 
substitutions  in  the  equation  follow  the  rule  that  sensors  a  and  c  must  be  at  right  angles  to  each 
other.  When  all  four  sensors  are  valid,  the  invariants  that  best  fit  the  four  sensors  readings  are 
given  by 


+  ab 


(7) 


(8) 

where  Ga,  Gb,  Gc,  and  Gd  represent  the  geophysical  stress  resultant  derived  from  four  sensors 
spaced  at  45°  intervals. 
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Figure  1.  Stress  components  at  Buoy  2  on  SIMI  Floe  4  expressed  in  lead  coordinates. 


Shear  invariant  (kN/m) 


Figure  3.  Relative  orientation  of  Buoy  3  and  other  stress  sensors  in  Lake  Andy  near  SIM  camp 
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Figure  4.  Stress  invariants  from  Buoy  3  (sensors  b,  c,  and  d)  (three-hour  averages  for  December  1993  to  February  1994). 


Shear 

invariant 

(kPa) 


Figure  5.  Stress  invariants  from  Buoy  3  (sensors  b  and  d)  and  Lake  Andy  site  (sensor  e)  (three-hour  averages  for  December  1993  to  Febmary 
1994). 
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Figure  6.  Stress  invariants  from  Buoy  3  (sensors  b  and  d)  and  Lake  Andy  site  (sensor  f)  (three-hour  averages  for  December  1993  to  February 
1994). 
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Figure  7.  Stress  invariants  from  Buoy  3  (sensors  b  and  d)  and  Lake  Andy  site  (sensor  g)  (three-hour  averages  for  December  1993  to  February 
1994). 
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Figure  11,  Indentation  pressure-area  curve  with  SIMI  stress  d 


